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ABSTRACT 

We present a unified model for molecular core formation and evolution, based 
on numerical simulations of converging, supersonic flows. Our model applies 
to star formation in GMCs dominated by large-scale turbulence, and contains 
four main stages: core building, core collapse, envelope infall, and late accre- 
tion. During the building stage, cores form out of dense, post-shock gas, and 
become increasingly centrally stratified as the mass grows over time. Even for 
highly-supersonic converging flows, the dense gas is subsonic, consistent with 
observations showing quiescent cores. When the shock radius defining the core 
boundary exceeds R ~ 4a(47rGpmean)~^''^5 where a is the isothermal sound speed, 
a wave of collapse propagates from the edge to the center. During the building 
and collapse stages, density profiles can be fit by Bonnor-Ebert profiles with tem- 
perature 1.2 - 2.9 times the true value, similar to many observed cores. As found 
previously for initially static equilibria, outside-in collapse leads to a Larson- 
Penston density profile p ~ 8.86a^/(47rGr^). The third stage, consisting of an 
inside-out wave of gravitational rarefaction leading to p oc r~^/^, v oc r~^/^, is also 
similar to that for initially-static spheres, as originally described by Shu. We find 
that the collapse and infall stages have comparable duration, ~ tff, consistent 
with estimates for observed prestellar and protostellar (Class O/I) cores. Core 
building takes longer, but does not produce high-contrast objects until shortly 
before collapse. The time to reach core collapse, and the core mass at collapse, 
decrease with increasing inflow Mach number. For all cases the accretion rate is 
^ a^/G early on but sharply drops off; the final system mass depends on the 
duration of late-stage accretion, set by large-scale conditions in a cloud. 



Subject headings: ISM: clouds — ISM: globules — stars: formation 
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Introduction 



Dense molecular cores are the immediate precursors of new stars on small scales, and 
understanding how they grow and evolv e is fundamental to the theory of star formation 
( IShu et al.lll987l : iMcKee &: Ostrikerl 120071 ) . Because many elements are involved in core for- 
mation, complete theoretical models have not yet been developed, and it is not yet clear 
whether a single dynamical effect dominates the overall process, or whether several con- 
tributing effects have comparable importance. In one limit that has been studied in some 
detail, ambient velocities are negligible, and self-gravitating cores form by the slow diffusion 
of partially-ionized gas thr ough strongly-supporting mag netic fields until a supercritical con- 
figuration is reached (e.g. iMouschovias fc Cioleklll999l ). In another limit, which has been 
considered more recently - but in much less detail at core scales, magnetic support is negligi- 
ble, and supersonic turbulence creates and destroys condensations, with some fraction of this 
gas s ufficiently dense and long-lived that it can undergo collapse (e.g. iMac Low fc Klessen 



20041 ). As observed clouds are both magnetized and strongly turbulent, the eventual theory 
for core formation that is deve loped must account for both processes; pio neering work to- 
wards this goal has begun (e.g. iKudoh fc Basull2008l : iNakamura fc Lill2008l ). Because of the 
technical challenges involved in building comprehensive models and the need to elucidate the 
contributing physics, it is important to develop simplified models in greater detail. In this 
contribution, we consider aspects of core growth and evolution in the turbulence-dominated, 
unmagnetized limit. 

Increasingly detailed observations in recent years provide constraints on theoretical mod- 



els (s ee e.g. the reviews of ldi Francesco et al.ll2007l : IWard- Thompson et al.ll2007l : lAndre et al. 



2008|). 



dimensional (angle-averaged) density profiles (e.g. IWard- 
2000l: iBacmann et~aDl200ol : ISvans et al.boOll 



fe.g. 


Ward-^ 


Thompson et al. 


1994; 


Shirley et al. 




Alves et al. 


2001 


Kirk et al 


boos 


Kandori et al. 



20051 ) generally show a uniform-density center surrounded by a power-law envelope extend- 



ing to an outer radius 



' Bonnor 


1956; 


Ebert 


1955) 



The interpretation in terms of static equilibrium is problematic, however, insofar as many 
cases show center-to-edge density contrasts exceeding the maximum ratio {pd pedge = 14.0; 
here pc is the central density) that would be stable against collapse, and would also re- 
quire central temperatures greater than observed values in order to provide support for the 
total masses inferred from the integrated continuum emission. In addition, cores are gener- 
ally not isolated; rather than being surrounded by a high-temperature, low-density medium 
with pressure matching the core's outer edge, they are surrounded by mo derate- density 



cold molecular gas representing clumps and filaments within larger clouds (IBergin fc Tafalla 



20071 ). The interpretation of observed density profiles as static solutions is also not unique. 
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in that dynamically-collapsing cores initiated from near-c ritical equilibrium show the same 



density profiles as (supercritical) static solutions (see e.g. lMyersll2005l : iKandori et al.l 12005 
and below) . Concentra tions formed within turbulent flows can also have density profiles 
resembling BE spheres ( iBallesteros-Paredes et al.ll2003l ). 



Velocity information can distinguish between static, oscillating (IKeto et al. iBroderick et al 



20071 ). and collapsing cores, and can potentially also help discriminate how these cores formed 
out of more diffuse gas. Dense, low-mass cores generally have subsonic inte rnal velocity 
dispersions, whether f or isolated cores o r for cores found in clusters ( e.g. 



Goodman et al. 1998: Caselh et al. 2002; Kirk et al. 2007: Andre et al. 2007; iLada et al. 



MversI Il983l: 



20081 ). In cores containing protostars, signatures of infall on small scales (~ 0.01 — 0.1 pc). 



of molecular lines that trace high-density gas 



[e.g. 



Zhou et al. 


1993: 


Mardones et al. 


1997 



Gregersen et al.lll997l : Idi Francesco et al.ll200ll ). For starless cores, inward motions are often 
evident over both small ( ^ 0.1 pc) and larger scales 0.1 — 0.4 pc), sornetime s encom- 

see also 



passing a whole st ar-forming complex (ILee et al.l 1 19991 . l200ll : IWalsh et a. 



Peretto et al.ll2006l ). Small-scale inward motions within co res are subsonic ( 



2006 



Lee et al.ll200ll ) 



while larger-scale motions can be transonic or supersonic JWalsh et allbooef infer higher ve 



locities in lower-density gas), and may be indicative of converging larger-scale flows in which 
dense gas builds up in a shock-bounded stagnation region. 

The relative durations of starless and protostellar (i.e. containing an accreting embedded 
Class or I object) stages of core evolution are determined by comparing the relative numbers 
of the two classes of sources in a given cloud. Absolute core lifetimes are further obtained 
by comparison to the number of T Tauri stars with measured ages. Several studies using 
this statistical approach in different clouds have reached similar conclusions: the durations of 



the starless and accreting stages of cores are comparable (IBeichman et a. 



1999l : ljessop fc Ward^rhompsonlboool. 



200 



i 



Enoch et al. 



2008 



Evans et al 



Kir 



c et al.ll2005l : I J0rgensen et al. 



1986 



Lee fc Mvers 



2007; Hatchell et al. 



l2008h . Typical starless core lifetimes are estimated at 



2 — 5 X 10^ yr, amounting to a few times the free-fall time 



\32GpJ 



1/2 



1.37 X 10^ yr 



10^ cm 



-1/2 



measured at the mean core density p = lAmHUH- With lifetimes considerably below the 



ambipolar diffusion time for strong magnetic fields Iad ~ lOt// (e.g. iMouschovias fc Ciolek 
19991 ). this suggests that observed cores are gravitationally supercritical with respect to the 
magnetic field. This conclusion is also supported by magnetic field Zeeman observations, 
i ndicating that cores have mean mass-to-magnetic-flux ratios two times the critical value 
( iTroland &: Crutcherl l2008l ). Since cores are only identified in millimeter and submillime- 
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ter continuum when the uh exceeds a few xlO'^cm ^, in principle it is possible that an 



lence accelerates amb i polar diffusion, how ever (e.g 



Zweibel 


2002: 


Fatuzzo & Adams 


2002 



Li &: Nakamural |2004J : iHeitsch et al.l |2004| ) , so it is also possible that the flux loss needed 
to reach a magnetically supercritical state may occur more dynamically, at densities below 
10"^ cm~^. 

Theoretical modeling of core evolution has a long history. Much work has focused 
on the evolution of unstable thermally-supported equilibria into collapse (formally result- 
ing in infinite density at the origin), followed by accretion of the envelope. Self-similar 



solutio n s for col l apse a n d /or accretion sta ges of i s otherm al spheres were found by iLarson 
/Il969h. IPenstonI ( 196£ [). Shul (119771). and iHunterl (|l977n: these \ vere later generalized by 
Whitworth &: Summers ( 1985 ). Larson ( 19691 ) and Penston ( 1969 ) (hereafter LP) indepen- 
dently found self-similar solutions which describe the density and velocity prior to the instant 
of protostar formation (defined as the instant at which the central density become infinite). 
In the LP solution, the radial velocity approaches a constant value —3.28a and the density 
approaches 



p = 8.86- 



at the instant of central protostar formation, with mass inflow rate M = 29.1 a^/G. Here, a 
is the isothermal sound speed, and the dimensional factor in the accretion rate is given by 



^ = 1.6 X 10"^ M0 yr- 



/ rp \ 3/2 

V 10^/ 



(3) 



The analysis of IShul (119771 ) showed that for an initial profile that is a static singular isother- 
mal sphere, p = 2a^/(47rGr^), evolution yields an "inside-out" infall solution in which a 
wave of rarefaction propagates outward at the sound speed. Inside of the expansion wave, 
the mass inflow rate is M = 0.9 7 5a^ /G independent of r, and gas accelerates to free-fall 



[V oc r' 



-1/2 



p (X r 



-3/2^ _ 



Hunter! (119771 ) connected and extended the investigations of LP 
(which address evolution prior to protostar formation) with that of Shu (which focuses on 
the accretion stage). He showed that self-similar solutions before and after the point of 
singularity formation (i.e. t = 0) can be smoothly matched. This allowed the LP sol ution to 
be ext ended into the accretion phase with similar free-fall behavior near the origin; iHuntei 



(1977) also found a sequence of self-similar solutions valid for all time that approach the 
( 119771 ) expansion wave solution. 



Shu 



Many numerical simulations of isothermal collapse have shown that the density in the 
core approaches a p oc profile at the point of protostar formation, regardless of how col- 



lapse is initiated (IBodenheimer fc Sweigartlll968l : lLarsonlll969l : |Penstonlll969l : lHunterlll977 
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2003 



Vorobyov fc Basu 



2005 



Foster fc Chevalierl Il993l: lOgino et al.l Il999l: iHennebelle et al.l l2003l : iMotoyama fc Yoshida 



Gomez et al.l 120071 ). These simulations include initiation from 



a static configuration that is unstable, and initiation from static, stable configurations that 
are subjected to transient compression, either from enhanced external pressure or a con- 
verging velocity field. Another feature common to the results from simulations is that the 
collapse generally begins on the outside, with the infalling region propagating inward as the 
central density increases. At the time of singularity formation, the central velocity has been 
found to be comparable to the value —3.3a derived by LP, with the densi ty normalization 
also similar to the LP result (iFoster &: Chevalierl Il993l : lOgino et al.lll999l ). Following the 
instant of protostar formation, the evolution of the mass accretion rate over time depends 
strongly on the initial conditions, however. 



Simulations with triggered core collapse (IHennebelle et al.l l2003l : iMotoyama fc Yoshida 



20031 ) are motivated by the fact that star-forming regions are highly dynamic, such that ex- 
ternal compression may significantly affect core internal evolution, and enhance the accretion 
rate by raising the central density. Triggering events may be associated with high mass star 
formation, but even without these highly energetic events, the large-scale turbulence that 
pervades giant molecular clouds (GMCs) can compress initially-quiescent cores. Taking this 
idea one step further, it is interesting to consider not just th e core collapse proce ss, but also 
the core formation process, in a strongly turbulent medium. iGomez et al.l (120071 ) conducted 
one such study, considering how an impulsive converging velocity field can create gravitation- 
ally bound, centrally-concentrated cores. Core formation induced by supersonic turbulence 
has also been studied in a number of numerical simulations that focus on the large scales, 
with much of the emphasis on determining the distribution of core masse s for comparison to 
observ ed core mass functions and the stellar IMF (see e.g. the review of iMcKee fc Ostriker 
(120071 )). However, these studies have not had sufficient resolution to investigate the internal 
properties of the cores that form. If the mass of a core is built up over time as the post-shock 
product of colliding supersonic flows, what is the detailed evolution leading up to collapse, 
and during the accretion phase? 

In this paper, we initiate a study of dynamically-induced core formation and collapse 
by considering perhaps the simplest possible situation: a supersonic, converging, spherical 
flow. Our initial conditions are a uniform low-density medium with no stratiflcation. A 
dense core forms inside a spherical shock, and over time becomes stratifled as its mass grows 
and it becomes self-gravitating. When the stratiflcation becomes too great, collapse and 
subsequent accretion occurs in a similar manner to the case of an initially-unstable static 
equilibrium. We consider cases of varying inflow Mach number, and with the large-scale 
inflow either steady over all time, or shut off after an interval. 
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The condition of spherical inflow that we adopt for this first study is, of course, hkely 
to be rare in real clouds. As the main purpose of this study is to take the first step to- 
wards unified models of core formation and collapse in dynamic environments, however, we 
consider one-dimensional solutions the natural place to start. We shall show that many fea- 
tures consistent with observed cores are evident even in these idealized models, suggesting 
that they are generic to dynamic core formation scenarios. The present set of simulations, 
in addition to enabling identification of characteristic evolutionary stages, also serve as a 
useful reference point for more realistic but more complex simulations. More typical than 
a converging spherical fiow would be a converging planar fiow, which yields a dense post- 
shock stagnation layer in which self-gravitating cores can form. Numerical studies that we 
have begun for supersonic planar infiows show results for core building and collapse that are 
qualitative similar to the present results for supersonic spherical infiows. 

The plan of this paper is as follows: In §2 we present the governing equations and 
describe our problem specification and numerical method. Section 3 briefiy describes results 
of collapse initiated from static configurations, demonstrating that we reproduce prior results, 
and providing a baseline for comparison to our models of dynamic formation and collapse. 
Section 4 presents numerical results for our converging-fiow simulations, covering the stage 
of core formation and evolution up to the point of singularity formation in §4.1 and the 
subsequent stages in stages in §4.2. We introduce a breakdown into new physically-defined 
stages in §4.3, and quantify the evolution of accretion rates in §4.4. Section 5 summarizes 
our new results and discusses our findings in the context of previous theory and observations. 



2. Governing Equations and Numerical Methods 



The equations of motion for a spherically symmetric fiow take the form: 

dp 1 d{r'^pv) 



dt 

dv dv 
dt Or 



dr 

IdP 
p dr 



GM{r) 



(4) 



(5) 



where M is the mass within radius r defined by dM = Anr'^pdr, v is the radial velocity, P 
is the gas pressure, and p is the density. For prestellar collapse, an isothermal equation of 
state P = a?p is often used, because cooling is so effic ient that the gas remains at n early 
constant temperature during the gravitational collapse (jLarsonlll969l : lNakamuralll998l ). We 
adopt an isothermal equation of state. 



For ease of comparison with previous work, we introduce the following dimensionless 
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variables: 

r = t/to, (6) 

^ = r/ro, (7) 

D ^ p/po, (8) 

u = v/a, (9) 

m = M/Mq. (10) 

Here po is a fiducial density representing the volume-averaged ambient density in a cloud on 

large scales, which we shall later set to the uniform density of the converging flow. The unit 
of velocity is the isothermal sound speed 

/ rp N 1/2 

a = 0.19 km s-^ (^^j , (11) 
the unit of time is related to the free-fall time at the fiducial density by 



the unit of length is related to the Jeans length at the fiducial density Lj = a(7r/G'po)^ by 



^° " (4.Gpo)^/^ = 2^^ = IiOkJ llOW^J ' ^^'^ 

and the corresponding basic unit of mass is po'^o = a^/[47r(47rG''^po)^^^]- The mass unit 
adopted in equation (fTOi) is larger than this by a factor 47r: 

«^ f T V'^ ( uh \ -1/2 



With the dimensionless variables, the fluid equations become 

du du 1 dD m 

m = [ Dedi. (17) 



W e solve the ID hydrodynamic equations f|T5|) - f|T7|l with the ZEUS-2D (jStone &: Norman 
19921 ) code, in spherical symmetry. During the evolution to singularity formation (collapse 
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phase), we adopt an inner reflec ting boundary condition. For the pos t-collapse accretion 
phase, we implement a sink cell ( iBoss &: Blacklll982l : lOgino et aljll999h at the origin when 
the central density reaches a reference value. Subsequently, the inner boundary condition 
is changed to outflow, and the value of the central point mass is tracked via the integrated 
flow off the grid, with M^r = [a^ / G)DinUinifn- The sink cell is only implemented after the 
inflow in the central region become supersonic, so that information from the inner boundary 
cannot propagate into the remainder of the grid. 



3. Evolution of Initially-Static Cores 

As discussed in §1, many previous numerical si mulations of core evolution have adopted 



static initial cond i tions and a flxed total m ass (e.g.. iHunterl 119771 : iFoster fc Chevalierlll993 



Ogino et al.l Il999l : IVorobyov fc Basul l2005l ). The adopted initial density proflle shapes are 
consistent with (or similar to) a hydrostatic equilibrium, i.e. a BE sphere, with the ini- 
tial density perturbed above the value that can be supported by the internal pressure in 
order to initiate collapse. These flxed-mass simulations adopt a pre scribed extern al pres- 
sure at a low density with an effectively Lagrangian outer boundary (IHunterl Il977l). or else 



a flxed outer boundary with p rescribed external p r essure and no mass inflow (lOgino et al. 



19991 : IVorobyov fc Basul l2005l ) . iFoster fc Chevalierl (119931 ) explored both types of boundary 
conditions, and found very similar results for a given initial cloud density proflle. They 
concluded that the evolution for the flxed-mass case is insensitive to the outer bound- 
ary condition for initiall y -unst able equilibria. This is consistent with the argument of 
Bodenheimer fc SweigartI (119681 ) that the outer boundary condition does not affect evolu- 
tion up to collapse as long as the free-fall time (eq. HD is shorter than the cloud crossing 
time. The ratio of the free-fall time at the mean density to the sound crossing time r^ax/O' 
over the radius of a BE sphere is vr/ (8^^)^^^, which approaches 7r/4 ~ 0.785 as ^ approaches 
infinity. Here, \E' = the dimensionless gravitational potential. For the critical case, the 

free-fall time is 0.71 times the sound crossing time over the radius. 

For comparison to previous work, we consider collapse of an initially-static BE sphere. 
For initial conditions, we adopt a critical BE sphere, i.e. the outer boundary of the grid is 
at radius 



R 



BE,crit 



6.45 



(47rGp,)i/2 



1.72 



{AnGp 



edge 



1 1/2 



2.70 



(47rG'p„)^/^ 



;i8) 



corresponding to dimensionless outer radius = ^cru = 6.45(po/Pc)^^^- Here, Pmean is the 
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total core mass divided by its volume. The mass of the critical BE sphere is 



BE,crit 



1.18 



4.18 



(47rG3p,,,,)i/2 



1.5 Mp 







T 



10 K 



3/2 



^H,edge 

W cm-^ 



-1/2 



(19) 

To initiate collapse, density is perturbed above the equilibrium value by 10%. Our outer 
boundary condition is at a fixed pressure, with no infiow. The temporal evolution of the 
accretion rate for this model is shown in Figure [TJ 



Features similar to those outlined by IVorobyov &: Basul (120051 ) are observed in our sim- 
ulation. The accretion rate peaks at a value approaching the LP prediction at the in- 
stant when the central d ensity (formally) becomes infinite, and steeply decreases thereafter. 
Vorobyov fc Basul (120051 ) find that the decline in the accretion rate after singularity forma- 
tion (unlike the increase in the accretion rate in Hunter's self-similar extension of the LP 
solution) can be attributed to the variation of velocity with radius in the sphere as it evolves 
toward collapse. We also considered cases of much larger initial static spheres, with outer 
boundaries = 5 ^crit corresponding to highly-unstable BE configurations. For these cases, 
the accretion rate decreases until it reaches a plateau at M = 1.45 0^/(7, consistent with the 
value reported in IShul (119771 ) (when the density is 10% greater than for hydrostatic equihb- 
rium), and then further declines to zero after a rarefaction wave propagates inward from the 
outer boundary to reach the center. 

Similar to the results of previous simulations, we find that the first collapse is "outside- 
in," with velocities initially nonzero only in the outer parts where the imbalance between 
gravity and pressure is largest^ This is because the inner portions of the sphere, at ^ < icritj 
initially are equivalent to stable BE solutions. In all of our models initiated from static 
spheres, the density profile approaches the LP self-similar solution D = 8.86^~^, and the 
velocity in the inner region approaches —3.28a, at the moment of singularity formation. 
Before this time, the density profile in the central region is fiat with a magnitude that 
increases over time. The process can be thought of as a wave of compression propagating 
from the outside to the inside, creating a density distribution in which the ratio of radius to 
the Jeans length at the local density is everywhere constant: r/Lj(r) ^ V8.86/27r ^ 0.47. 
The singularity represents the instant the compression wave converges at the center. 

Supersonic infiow velocities can be achieved without shock formation in the interior of 
the core (except at r = 0) because inward acceleration occurs at all radii where the in- 
ward gravitational force exceeds the outward pressure force. By construction, the solutions 



^For cases where the init ial sphere is larger than th e critical BE sphere, collapse begins at radii near £,crit, 
as was previously shown by [Foster fc Chevalieii (jl993( ). 
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initiated from static configurations all have gravity (slightly) exceeding pressure forces every- 
where in the initial state, so that the inward acceleration is nonzero. After a singularity forms 
at the center, accretion begins, and the flow in the interior transitions from the v = const., 
p oc LP solution to a free-fall solution, with the accreting reg ion propag ating from the 
interior to the exterior in a manner similar to that described by IShul (119771 ). Thus, while 
collapse develops in an "outside-in" fashion, accretion develops in an "inside-out" fashion. 
Note that accretion in any centrally-concentrated configuration should work its way out- 
ward from the center, because gravitational collapse times decline outwards oc p^^^^, which 
is t// oc r for an inverse- square density profile (true for either the LP profile or the singular 
isothermal sphere). 



4. Converging-Flow Model Results 



In this section, we present the results of our simulations of core formation and evolution 
for the case of converging spherical, supersonic flows. For these simulations, the outer 
boundary condition at ^max is maintained at a constant density and inflow velocity. The 
inflow velocity is characterized by the Mach number relative to the isothermal sound speed. 
Ma = Vin/a. We vary over the range from 1.05 to 7. For some models (see below), we 
suppress inflow at the outer radius after collapse occurs to test how the late-time accretion 
rate is altered. The initial conditions consist of uniform (low) density, and uniform inflow 
velocity equal to the value at the outer boundary. The size of the grid, in terms of the 
reference length scale given in equation (fT3!) . is 2.51327ro, which amounts to a fraction 0.4 
of the Jeans length at the initial density. Note that the radius of a critical BE sphere at the 
same external pressure (so that pedge = Po) would be RBE,crit = 6.45(ped£,e/Pc)^^^^o = l-72ro. 
That is, a critical BE sphere confined by the same ambient pressure would be able to fit 
within our simulation domain, with plenty of room to spare. The size of zones in the radial 
direction has a constant logarithmic increment, i.e., Arj+i = (1 -|- a)Ari, for some a > 0, 
such that Ti = ri + a~^[{l + a)*~^ — l]Ari and Ari/(rjnax — '"mm) = "[(1 + a)^~^ — 1]""*^. 
For all the converging flow simulations in this paper, a is set to 0.009; 605 and 597 grids are 
used during the collapse phase and accretion phase respectively. 



4.1. Core Formation and Collapse 

We begin with a description of the core formation process, which is similar for all of the 
converging-flow models. Because of the reflecting boundary condition at the center and the 
initial inflow velocities, immediately after we initiate the simulation, a shock forms at the 
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origin and propagates outward. The inflowing matter is compressed by the shock. The shock 
front divides the converging inflow into two regions: an inner dense post-shock region and an 
outer low-density region of supersonic inflow. These two regions evolve quasi-independently 
but are connected by shock jump conditions. Under competition between gas pressure and 
self-gravity, the inner region contracts slowly to begin forming a dense core. As self-gravity 
starts to overwhelm gas pressure, the dense core enters the collapse phase. 

At the beginning, the inner region is quasi-hydrostatic, with the velocity a linear function 
of radius. With negligible density gradient and self-gravity in the early stages, equation (HM 
becomes 

^+.^«0. (20) 

This equation is satisfied by 

u = (21) 
r - To 

for \t/tq\ ^ 1, m ^ ~^/to; i-e. the coefficient of the linear profile is constant in time. The 
leftmost lower panels of Figures [2] -[5] show this linear- velocity behavior in the shocked region, 
at the time that equals half of the collapse time Tcou- Throughout this paper, we define the 
collapse time Tcou as the moment that pd po = 4 x 10''. Shown in the leftmost upper panels 
of Figures [2] - [5] are the density profiles; even when the density profile has nonzero gradients, 
the velocity in the inner region is dominated by the linear term. Notice that the leftmost 
panels have linear-linear scales. The velocity over the whole post-shock region is subsonic 
and negative, which means the core is slowly contracting. For the = 4, 7 models (see 
Fig. m [5]), the inner part of the density and velocity profiles oscillate at the beginning of 
simulations. As the shock front propagates outward, the mass inside the shock increases, 
and so does self-gravity. After a period of accumulation lasting about 90% of the time until 
collapse Tcoii, the slowly-contracting dense region starts to be gravitationally unstable. 

Note that the density and velocity profiles of the regions outside of the shock go through 
a transient evolution after simulations begin. The density profile evolves from a uniform 
profile set in the initial condition to a p oc profile consistent with supersonic radial inflo\\EI; 
the material is also increasingly accelerated by gravity. In Figures [2] -0, the inflection feature 
in the velocity profiles exterior to the shock corresponds to a wave propagating inward at a 
speed equal to the inflow speed plus the sound speed. 

Since post-shock velocities are subsonic, gravitational instability is expected to develop 
at a point when the radius of the inner dense post-shock region becomes comparable to the 



^For supersonic radial flow, v ~ const, upstream of the shock combines with the steady state mass inflow 
condition pvr^ — const, to yield p cx r^^. 
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critical radius of a BE sphere. This expectation is indeed borne out by our simulations, which 
moreover show that the properties of the collapsing dense inner region are similar as those 
of collapsing cores initialized from hydrostatic BE spheres. The collapse follows an "outside- 
in" pattern, starting from the shock front. The central density increases dramatically and 



Larson 


1969; 


Penston 


1969; 


Hunter 



19771 : IShulll977l ). which approaches the LP density profile D{^) = 8.86^ ^, and the velocity 



approaches —3.28a. 

The three panels on the right of Figures [2] - [5] show the density profile (top) and the 
velocity profile (bottom) at three different instants during the outside-in collapse. The 
density profiles plotted are normalized by the central density and the interior regions are fit 
by BE sphere density profiles. % The ratio between the fitted BE sphere temperature Tbe 
and the true temperature Tq is also noted in the figures. The first figure of these three shows 
profiles at the instant when the radius of the post-shock region reaches the critical radius (see 
eq. [18]) of a BE sphere at temperature Tq with the same central density i.e. rghock = RBE,crit- 
We define the period after this as the supercritical regime. The center-to-edge density ratio 
after this point exceeds ~ 10. The second figure of these three shows the instant when the 
radius of the post-shock region reaches twice the critical radius of a BE sphere. The third 
figure is at the instant the core collapse Tcou (defined here as Pc/po = 4 x 10''). The long 
dashed diagonal line in the third figure shows the LP density profile, which is very close 
to the numerical solution. The time interval between r shock = RBE,crit and collapse for the 
Ma = 2 case in code units is 0.048 (see eq. [T^ . which corresponds to 1.1 x 10^ yr if the 
inflowing ambient medium's density is uh = 100 cm^^. 

We note that prior to collapse, the velocities in the dense gas (inside the shock) remain 
small. In particular, for Ma > 2 cases, the inflow velocity inside the shock is subsonic 
throughout the post-shock region before the time when rghock = RBE,crit- When Ma is high. 
Figures H] and [5] show that the post-shock velocities remain subsonic up to the instant of 
collapse. 

The outward propagation speed of the shock, which from the simulations is approxi- 
mately constant (Fig. at early time, can be obtained using the shock jump conditions. If 
the shock position is ^sh = UshT, then provided equation (pTj) is satisfled, the velocity on the 
downstream side is 

u, = (22) 

T -To 



■^That is, we match solutions of the density profile inside the shock to solutions of the hydrostatic equation 
= —^^^1^ with a fitted temperature Tbe- These fits are further discussed below. 
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Notice that in order for the downstream velocity to be inward and the shock to be propagating 
away from the origin, we must have r < tq. The isothermal shock jump conditions, with 
subscript "m" denoting upstream and "rf" downstream values, aie Dd{ud—Ush) = Du{uu—Ush) 
and Ddll + iud-UshY] = -D„[l + (m„ -Ms/^)^]; together these imply Dd/D^ = {uu-Ushf and 
{ud—Ush){uu—Ush) = 1- Treating the shock speed as approximately constant so that equation 
(!22|) holds, we can solve for the shock velocity to obtain Ush = ^{uu + + 4 — 4(r/ro)]^''^). 
For t/tq ^ 1, and taking the upstream velocity as m„ ~ —Ma which corresponds to the limit 
of strongly supersonic inflow, gives 

u^f^ ^ - \-Ma + iMj, + 4)1/2] _ ^ ^23) 

2 L « ^ V a ^ ; J Ma + {Ml + 4)1/2 V ; 

The shock speed decreases as the Mach number increases, and therefore from equation fl22l) 
the post-shock flow speed also decreases as the Mach number increases. Figure [6] shows the 
position of shock front versus time for Ma = 4. The intercept and the slope are based on the 
best fit of the linear part where r G [0.0, 0.04]. The intercept is nearly and the slope 0.3061 
is the measured shock speed in units of the isothermal sound speed a; equation (123|) predicts 
a slightly smaller value 0.24. The analytical solution (eq. 123!) for function of Ma is 

plotted as a dotted line in Figure [7] (labeled as Ush,estimate) and the shock speeds measured 
directly from simulations Ugh are plotted as asterisks. The analytical approximation is about 
15-28% below the measured value from the simulations as Ma ranges from 1.05-7.0. 

Using the constant-shock-speed approximation, the immediate post-shock density can 
also be obtained in terms of the pre-shock density, in the limit t/tq ^ 1 and using u„ ^ —Ma 
as: 

^^ = \[Ma + {Ml + AYl^]\ (24) 

which for Ma ^ 1 is Dd/Du ~ M^. Because of the radial convergence of the inflow in 
the simulations, varies; it is initially equal to 1, but after an initial transient, in the 
highly-supersonic limit Du would approach {pouter / ^shY because of mass conservation. 

It is interesting to investigate how the state of the core when it collapses depends on 
Mach number. Figure [H] shows the density profiles and the velocity profiles of the simulated 
cores at the instant of collapse Tcoiu for different Mach numbers. Firstly, it is clear that all 
the density profiles approach the LP solution, especially for low Mach number cases. The 
inflow velocity inside the shock is supersonic and does not strongly vary over the dense core 
region. The smaller the Mach number is, the closer the inflow speed is to the LP result, 
V = —3.3a. Secondly, the radius of the core at the instant of collapse decreases as Mach 
number increases. We plot this dependence in Figure [7] with diamonds. 

We quantitatively compare the basic core properties in Figure U\ which shows the col- 
lapse time Tcoii and the shock radius C,coU at time Tcou, both as function of M^. As noted 
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above, ^coii decreases with increasing M^, the same is true for Tcou- We also recall that Ush 
decreases with M^. Interestingly, while ^coihUgh, and Tcou all decrease with Ma, the ratio 
icoul iush'Tcou) is nearly constant with M^: it ranges only from 0.34 to 0.42 (see Fig. [7]). 
This result is potentially useful for empirical estimates of core lifetimes, since the inflow 
velocity Mq, the isothermal sound speed a and the radius of dense core are all in principle 
measurable. If the ratio icoiil {ushTcoii) is taken as a constant ~ 0.4, and the shock speed is 
estimated via equation (123|) . then the life time of cores up to the point of collapse is given in 
dimensional form by tcou ~ ^■'^RcoiA^a + + 4)^/^]/a. In practice, it may be difficult to 
measure Ma outside of a core, because the density is much lower than that of the core, and it 
is difficult to isolate the immediate environment of the core from foreground and background 
gas. This result is still useful in a statistical sense, however, using the mean Mach number 
of the turbulent flow in a cloud. 

Another direct observable is the core density, so it is interesting to test how the values of 
the collapse time and radius depend on the mean density in the core at the time of collapse. 
Diamonds in the top panel of Figure [9] show the core radius in units of a{A'KGpmean)~^^'^ ■, 
which is C.coiiipmean/poY^'^, as a function of Ma. We can see ^coiiiPmean/ PoY^'^ is nearly 
constant, ranging from 4.58 to 3.42 as Mach number increases from 1.05 to 7. Taking this 
as approximately constant, and taking the measured core mean density, the predicted size 
of core at the time of collapse is Rcoii ~ 4a(47rGpmean)~^''^- Note that this radius is ~ 50% 
larger than the critial BE radius for the same temperature (see eq. [TSl) . Since the post-shock 
density increases relative to the upstream density approximately as Dd/Du ~ Ml (see eq. 
UM . it is also interesting to test how icouMa depends on Mach number. In fact this quantity 
decreases with Ma, as seen in Figure [E 

To express the core collapse time in terms of observables, we normalize the collapse time 
using the mean core density. This quantity tcoZ«(4vrG'p„ean)^^^ = Tcon{pmean/ PoY^'^ , is plotted 
in Figure [HI as a function of M^. For reference, the core collapse time normalized using Ma is 
also plotted in Fig. [91 The free-fall time t// for a uniform sphere in units of (47rGpmean)~^^^ 
is 1.92, so that we have tcoii ~ 8 — 26 t// as Mach number varies from 1.05 to 7. This time 
scale is much longer than the observed values ~ 2 — 5 t/j for prestellar cores. The reason 
for this disparity is that during the early part of its evolution, the central density of the 
core is low, and it would not be identiflable within its surroundings. This is evident in the 
low contrast between the center and the edge of the core seen in the flrst frames of Fig. [2] 
- [5l Only when the core approaches collapse does the center-to-edge contrast become large. 
For example, for the Ma = 1.05 model, the center-to-edge density ratio pd Pedge reaches 2 at 
r = 0.78 and reaches 5 at r = 0.92; this can be compared to the total time until collapse, 
TcoU = 1.027 for this model. For the whole set of models, the observable fraction of the pre- 
collapse core life is 30 — 50% if we choose pc/ pedge > 2, or 10 — 20% if we choose pd pedge > 5. 
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Taking the period when pd Pedge > 2 or pd Pedge > 5 as the period over which a core could 
be observable, e.g. in submilimeter continuum, our simulations give tobservaUe ~ 2 — 13t// or 
tobservabie ~ 1 — 6t//, respectively. The latter is consistent with observed estimates. 

Another observable aspect of prestellar cores is their density structure. As discussed 
above, for each solution at the times shown in Fig. [2] - [5l we fit a BE sphere profile. The 
density in code units po is -0(0 = p/Po and the radius in code units is ^ = r^A-nGpQyf^ / a. 
When the central density is instead used to normalize, the density and radius variables are: 

D = p/p, = Dpo/pc (25) 

and 

i = r = UPc/Po ' ■ 26 

a 

For a BE sphere with sound speed a be, the density profile normalized by the central density 
is Dbe, which is a function of the scaled radius ^be = r{4:TrGpcY^'^ /ube = ^o./ c^be- To fit 
the density profiles in our simulations to a BE sphere, the only free parameter is the ratio 
of the isothermal sound speed a^E of the BE sphere to the sound speed in the simulations, 
a. Thus, for any given density profile D we adjust the value of aBE/o. until a good match 
between Dbe and D is obtained. This procedure yields the fitted temperature ratio: 

^ = i^)" . (27) 

The values obtained for Tbe/Tq are marked in Figures [2] - [5l The range of values we find 
is 1.23 to 2.89. This range is consistent with theoretical expectations. As the radius of 
the BE sphere extends to infinity, the density profile approaches the singular solution p = 



2a|^/(47rG'r^) (IShul 119771 ) corresponding to Dbe = p/Pc = 2{aBE/aY^~'^- The density 
profile of the collapsed core approaches the LP solution D = 8.86^~^. To match the LP 
profile with the singular profile therefore requires 2 {aBE/o)'^ — 8-86, which corresponds to a 
temperature ratio Tbe/Tq = 8.86/2 = 4.43. Fits of core profiles to BE spheres that indicate 
values of Tbe up to 4.43 times the measured thermal temperature therefore are expected if 
collapse has taken place; this does not by itself indicate that magnetic support is present. 



4.2. Post-Collapse Evolution: Infall and Accretion Stages 

After the central density becomes singular, the evolution transitions to the infall and 
accretion stages. For our simulations, we make this transition by implementing an outflow 
boundary condition at the center when the central density reaches 4 x 10 ''po- The initial mass 
of the central protostar is calculated by integrating the innermost part of the density profile 
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where density is between [1 x 10^,4 x 10^] po- The specific choice of this density does not 
significantly affect Tcou, icoiu or the subsequent evolution since the central density increases 
dramatically only at the very end of the collapse stage. 

At the beginning of the accretion phase, the material inside the s hock falls ori to the 



protostar. The material approaching the protostar is in a free-fall state (IHunteii 119771 ) . The 



region of unsupported infall st arts f r om th e center and propagates outward, similar to the 



"expansion wave" described by IShul (jl977l ). The density profile inside the rarefaction wave 
changes from p oc to r~^/^ and the velocity profile changes from f oc to r~^/^. For our 
simulations, this infall stage ends as the infall rarefaction wave arrives at the shock front. 
This generally occurs very rapidly (in less than 10% of Tcoiu see below). 

For an initially static density profile Ar'"^, where A = Ka^ / [At^G) and i^' is a constant, 
if the gas pressure is negligible the average speed of the rarefaction wave is {2\/2K /t^) a. For 
K = 8.86, which is the LP profile, this yields 2.7a. For the real case, the initial velocity is 
nonzero and the gas pressure is non negligible, so that the rarefaction wave propagates at 
a modified speed. For example, for the Ma = 1.05 model, which has ^coii = 0.29 and infall 
interval Arj„j = 0.092, the average speed is 3.15a. For Ma = 4.0 and 7.0, the measured 
average infall speed of the infall rarefaction wave is 2.22a and 1.95a respectively. 

After the infall rarefaction wave arrives at the shock front, the final accretion stage 
begins, with material initially outside the shock falling onto the protostar at supersonic 
speeds. This process is similar to Bondi accretion, except that the central mass is growing and 
the velocity field for our simulations is uniformly converging at large distance. During this 
stage, the density and the velocity profiles vary oc r~^/^ and r~^/^ respectively, corresponding 
to free-fall. 

The typical density profiles and velocity profiles during the accretion stage for Ma = 1.05 
and 4 models are shown in Figure [TOl Three different points in the evolution are shown: 
the instant of core collapse, the instant when the infall rarefaction wave arrives at the shock 
front, and a point during the late accretion stage. The transition from the LP profiles to the 
free fall profiles in density and velocity are clearly evident in the figure. 

In a real system, the duration of the accretion stage depends on the environment of the 
protostar, and how long the inflow that creates the core is maintained at large scales. To 
explore how the late-time evolution is affected by changes in the accretion, we have conducted 
additional simulations in which the flow inward from the outer boundary is halted at the end 
of the infall stage (i.e. when the rarefaction reaches the shock). Suppression of inflow will 
affect the mass flux onto the protostar after the rarefaction wave from the boundary reaches 
the central protostar. We discuss comparison of these models to our standard models in 
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section 4.4. 



4.3. Definition of Evolutionary Stages 

Based on the results presented in §4.1 and §4.2, we have identified four main stages of 
protostellar core formation and evolution in a supersonic turbulent medium (see Figure [TT] 
for a schematic depiction): 

1. Core building — 

Converging fiows in a supersonic turbulent medium collide, with post-shock compressed 
gas accumulating over time in stagnant, shock-bounded regions. If these dense regions are 
not destroyed by larger scale turbulence, the high density gas will undergo a long contraction 
process during which gas pressure competes with self-gravity. The typical character of this 
stage is that the velocity inside the dense gas is subsonic and increases hnearly with distance 
from the center. Since the center-to-edge density contrast is relatively low, these clumps are 
gravitationally subcritical. Towards the end of this stage, when the center-to-edge density 
contrast becomes appreciable, these objects would become observable as prestellar cores. 
This core building stage lasts > 90% of t^oiu only the last 10 — 20% would be observable. 

2. Core collapse — 

The core built up during the first stage accumulates enough mass that it becomes grav- 
itationally supercritical, which we operationally define as r shock > RBE,crit- Self-gravity 
overwhelms the gas pressure, and the unstable core starts to collapse. The collapse begins 
near the shock front, where the imbalance between gravity and pressure gradient forces is 
greatest, and propagates inward. This collapse is an "outside- in" process. During core col- 
lapse, the central density increases dramatically and the infiow velocity becomes supersonic. 
As the collapse propagates inward, a density profile p oc is left behind and the velocity 
increases toward the center. The end of this stage is defined by the instant of protostar 
formation Tcoiu when the wave of collapse has reached the center of the core. The density 
profile approaches the LP profile, p = 8.86a^/(47rGr^). The velocity in the interior of the 
collapsed core approaches a uniform value comparable to that in the LP solution —3.28a. 
The central density is high enough to be optically thick and a protostar forms. Because 
of their large central volume densities, prestellar cores during this collapse stage would be 
observed as having high peak sub mm fiux densities. This stage lasts for a time Arsupcrit, less 
than 10% of the prestellar core lifetime Tcou- 

3. Envelope infall — 
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During this stage, the high density material inside the shock front falls to the protostar 
(or, if angular momentum were included, a circumstellar disk). This stage starts at the 
instant of the protostar formation and ends at the instant when the gravitational rarefaction 
wave reaches the shock front, clearing out the remnants of the dense gas that accumulated 
during core building. During this stage, the density and velocity profiles in the interior 
change from LP profiles to free-fall profiles. Since the core contains an embedded protostar, 
the system would be observationally classified as Class O/I. This stage lasts for a time Atj„j 
less than 10% of the prestellar core lifetime. In a real system, the fraction of the envelope 
mass that eventually reaches the center would depend on the details of the protostellar wind, 
which would sweep up at least the polar portion of the envelope, creating an outflow. 

4. Late accretion — 

During this stage, material from the ambient environment directly accretes to the proto- 
star (or, more realistically, a disk). Accretion during this stage is similar to Bondi accretion, 
with free-falling density and velocity profiles p oc r~^/^ and v oc r~^/^ over the whole core 
region. For a real system, the duration of the late accretion stage, and hence the final stellar 
mass, depends on ambient conditions far from the protostar. In a real system, the potential 
of this stage to contribute significantly to the stellar mass would also depend on protostellar 
winds, which can reverse the accretion. 

Although we have identified these stages based on idealized spherically sjnnmetric mod- 
els with constant gaseous inflow rates at large distance, we expect that the same stages 
would be present, in modified form, under more realistic conditions. Based on our simu- 
lations, the ratios Arsupcrit / TcoU and Arinf/Tcou decrease with increasing Ma, as shown in 
Figure [121 The supercritical stage and the infall stage have similar durations, and range 
from 9 — 3% of the prestellar core lifetime tcoii, which itself ranges from 8 — ^Gtjf. The su- 
percritical collapse period (stage 2) Atsupcrit and the infall period (stage 3) Atj„j thus both 
range over ~ 0.8 — Itff, as shown in Figure[T2l While the time to reach collapse would differ 
for nonspherical or nonsteady converging large-scale flows, we expect that the character of 
the evolution would not. We also expect that the ratios Atgupcrit/tff and Atinf/tff would 
remain order-unity. 



4.4. Evolution of Mass Accretion Rates 

Figure [13] shows the temporal evolution of the mass accretion rate and the total inte- 
grated mass of the central protostar for Ma = 1.05,2,4 and 7. These can be compared to 
the mass accretion rate and the integrated central mass for an initially static critical BE 
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sphere, as shown in Figure [T] (see also I Vorobyov &: Basul (120051 ) . who show similar accretion 



histories to Figured]). For both the initially-static collapse and our models that allow for 
core building from supersonic flows, there is a sharp early peak in the accretion rate. The rise 
starting from M = corresponds to the moment of protostar formation at Tcou- The smooth 
decline that follows (ending at the points marked "i" in Figure [T3l) corresponds to the infall 
stage, as the gravitational rarefaction wave propagates outward. At late times, however, the 
accretion differs for the initially-static vs. dynamically-built cores. For an initially-static 
unstable BE core, the late accretion steadily declines over time (Fig. [T]). In contrast, for 
cores formed in the Ma = 1.05 model (i.e. barely supersonic inflow), the late time accretion 
is nearly constant, and for cores formed in large Ma models, the late accretion rate increases 
over time (thin lines in Fig. [T3|) . 



The early-time peak accretion rates can be compared to the predictions of self-similar 
models. For the LP profiles at the instant of core collapse, D = 8.86^^^ a nd u = -3.28 



give M = 29.1a^/G, while the self-similar solution for the accretion phase in iHunterl (119771 ) 
predicts M = 46.195a^/G. In fact, we do see a jump in M above 29.1a'^/G as the evolution 
transits from the collapse stage to the infall stage. This phenomenon is most clearly evident 
for the Ma = 1.05 model, which has the highest resolution of the central region because the 
shock strength is lower than in the high Ma models, yielding a larger core (see Fig. [8]) at 
the instant of collapse. 

The detailed behavior of M during the late accretion stage can be understood in terms 
of various transitions that occur. For Ma = 1.05, the accretion rate (see in Fig. [T3]) starts to 
increase from point i until to point p, and then decreases. The increases from i to p occurs 
as gas stored between the shock front and the outer boundary collapses into the center. 
The point p represents the instant when the density in the whole outer region reaches a 
profile p oc r"'^/^. After point p, the gravitational rarefaction has reached the boundary, and 
subsequent accretion is limited by the inflow rate imposed at the outer boundary. For the 
Ma = 2, 4 and 7 models shown in Figure [131 the mass accretion process between point i and 
p is similar to that of Ma = 1.05. However, there is additional transitory behavior before 
the accretion rate decreases to the inflow rate imposed at the boundary. At the instant 
corresponding to point p, the rarefaction has produced p oc r~^/^ over the whole region. 
But the velocity profile has v oc r~^/^ only over the inner region. During the stage between 
point p and q, the density profile stays almost unchanged but the velocity profile evolves to 
reach v oc r~^/^ everywhere (see also Fig. [TOl) . After point q, the accretion rate decreases to 
the imposed inflow rate. The stage between p and q is most obvious for the Ma = 4.0, 7.0 
models. 



As mentioned above, we have also performed models in which inflow to the grid is halted 
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after the point when the rarefaction reaches the shock. The resulting late-stage accretion 
(see Fig. [T3l) is initially the same as in our standard models, but then declines over time, 
after the rarefaction wave reaches the boundary. 

The mass of the protostar at the end of the infall stage, and the total core mass Mcore 
inside the shock at instant of core collapse, are shown in Figure [TH For comparison, we also 
show the critical BE sphere mass (see eq. [19]) using the post-shock density at the time of 
collapse for Pedge] these are much lower than the actual core mass. Because there is continued 
mass passing through the shock during the infall stage, the post-infall protostellar mass is 
slightly higher than the core mass inside the shock at the time when core collapses. As 
the Mach number increases, the post-infall protostellar mass and the core mass at Tcou both 
decrease. The protostellar mass ranges over 0.06 — 8.8 Mq and the core mass at Tcou ranges 
over 0.05 — 7.5 Mq, taking uh = 100 cm~^ for the ambient density. 



5. Summary and Conclusions 

Star formation takes place in GMCs pervaded by supersonic turbulence, and theoretical 
models of prestellar (and protostellar) cores must take these large-scale supersonic flows into 
account. Here, we have developed models in an idealized, spherically-symmetric framework 
that nevertheless captures key aspects of the real situation, enabling us to identify and 
analyze the main stages of core formation and evolution in a dynamic environment. 

Our models differ from previous studies of core evolution in that the cores are not 
present as either stable or unstable density concentrations in our initial conditions - the 
initial density is everywhere uniform. Instead, during the first evolutionary stage cores are 
built "from scratch" by the collision of converging supersonic flows. The boundaries defining 
the outer edge of a core - where the density drops - correspond to a shock front across 
which the temperature is constant and the mass flux is nonzero. The shock front propagates 
outward, with the mass of the post-shock dense region growing in time. Initially, the core is 
essentially uniform. Over time, the mass grows sufficiently so that the core becomes centrally 
stratified due to self-gravity. Observationally, the latter part of this "core building" stage 
corresponds to prestellar cores that have low to intermediate peak brightness. The period 
over which Pc/ Pedge > 5, and a core would be clearly identifiable in observations, amounts to 
1 - 6 tff, with the free-fall time defined using the mean core density. 

When the center-to-edge density contrast exceeds ~ 10, the core becomes supercritical 
and a stage of violent "outside-in" collapse ensues. The density profile throughout the core 
approaches p oc r~^, and a protostar forms at the center. We define the instant that collapse 
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is complete and a protostar forms as tcoii, or Tcou in our dimensionless variables. Although 
the central density becomes very large, the wave of outside-in collapse still leaves most of 
the core mass in the outer parts. Observationally, this core collapse stage corresponds to 
prestellar cores that have high peak brightness. The period Atgupcru over which cores are 
supercritical, undergoing outside-in collapse, amounts to less than 10% of tcoii, or 0.8 - 1 tff. 

The third stage of evolution is governed by an "inside-out" wave of gravitationally- 
driven rarefaction propagating from the center of the core to the shock front that defines 
the core's outer edge. The accretion rate onto the star during this infall stage is initially 
very high, but declines over time. At the end of this infall stage, the dense envelope built 
during the first stage has plunged into the star. The velocity and density profiles approach 



free-fall, v (x r and p 



oc r 



-3/2 



respectively. Observationally, this stage corresponds to 



Class O/I embedded protostars. The period Atj„j over which cores undergo this inside-out 
infall is similar to the duration of the previous stage, Atsupcrit, and comparable to tff. 

During the final stage of evolution, there is no longer a massive envelope. The protostar 
can continue to accrete from the more distant, lower-density gas in its surroundings. The 
late-stage accretion rate and total mass accumulated by the system depend on the large-scale 
environment, rather than the properties established in the core during the building stage. 
Observationally, this stage corresponds to a non-embedded YSO that may still be accreting 
from a disk. 

Based on our simulations, our chief conclusions are as follows: 

1. The initiation of star formation via outside-in core collapse, followed by inside-out 
envelope infall, appears to be very robust. The dynamical behavior during these stages 
of evolution is very similar whether the core is initiated as an unstable equilibrium (as in 
previous models) or is built up dynamically through a shocked converging flow (as in the 
present work). The Larson-Penston singular solution with p = 8.86a^/(47rGr^) appears 
to be an "attractor," in that models initiated from stationary equilibria or with different 
supersonic converging velocities all arrive at this configuration at the moment of protostar 
formation. 

2. Prior to the point at which cores become supercritical and outside-in collapse begins, 
the velocities interior to cores are subsonic, even if they are created by highly supersonic 
flows. In fact, higher infiow velocities produce /ower post-shock velocities within the dense 
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3. Throughout both the core-building and core-collapse stage, density profiles for cores 
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formed by shocked converging flows can be fit by Bonnor-Ebert profiles, but with fitted tem- 
peratures Tbe larger than the true temperature T. The range of temperatures fitted for our 
models with Mach numbers up to 7 is Trr/T = 1 .2 — 2.9. This result is consistent with obser- 



vational findings (iKandori et al.ll2005l : iKirk et al.ll2005l ) that effective temperatures greater 
than directly-measured values (from fitting SEDs) are usually required in order to fit BE 
spheres to observed prestellar cores. The largest possible ratio that could be obtained for 
an isothermal spherical fiow is Tbe/T = 4.43, so that any observed ratio larger than this 
suggests that r nagne tic fields contribute appreciable support, or else the core is anisotropic. 
Dapp &: Basij (120091 ) have also recently pointed out that the temperature fit based on match- 
ing a BE profile may be significantly higher than the true kinetic temperature, for clouds in 
collapsing stages. 

4. At the time of collapse, for all Mach numbers the core size and mean density are 
closely related. We find that Rcou ~ Aa{ATrG pmean)^^^'^ within 15% for Mach numbers Ma = 
1.05 — 7. This radius is ~ 50% larger than the critical radius of a BE sphere with the same 
mean internal density. 

5. As Ma increases, and assuming a given ambient medium density po, the time to 
reach collapse tcoii is shorter, the physical size of the core at tcoii is smaller, the mean internal 
density at tcoii is higher, and the mass of the core at tcoii is lower. For high Mach number, 
the collapse time and collapse radius are related by tcoii ~ 2.6RcoiiMa/a. The range of 
core masses at the time of collapse at different Ma is consistent with observed core masses, 
although the specific dependence of core mass on Ma found in the present work may be 
sensitive to the spherical converging-fiow geometry we have adopted. 

6. The durations of the infall stage and collapse (supercritical) stage of evolution are 
comparable for all Ma, and are close to tff. This is consistent with observations indi- 
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(iBeichman et a. 



J0rgensen et al 
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2007 
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2008 : 



Evans et al. 20081). These 



stages are preceded by an extended core-building stage, during most of which the core would 
not be observable because its center-to-edge density contrast is low. 

7. The mass accretion rate onto the protostar (or, more realistically, star-disk system 
if angular momentum were included) peaks at the beginning of the infall stage at a value 
^ a^/G, and then declines steeply afterwards as the material stored in the envelope is ex- 



earlier work (see e.g.. 


Foster & Chevalier 


I993I: 
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Motovama & Yoshida 
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Vorobvov & Basu 


20051) for cores 



libria or which undergo externally-induced compression. Later accretion from the ambient 



-23- 



medium depends on how long the large-scale cloud maintains a focused converging flow. 

As noted above, some of our specific conclusions are hkely to change for non-spherical 
geometry, and for time-dependent rather than steady large-scale infiow. Furthermore, other 
elements that are present in real star formation have been entirely omitted in these models; 
these elements include rotation, which would lead to disk formation; protostellar winds, 
which would sweep up and remove a portion of the envelope during the infall stage, and 
could prevent late accretion altogether; and magnetic fields, which would alter the timescales 
and details of the evolutionary stages. 

We expect, however, that many of our basic results will carry over even if the ide- 
alizations we have adopted are relaxed. While large-scale supersonic converging flows in 
real GMCs are not generally spherical, the association of observed cores with high density 
surroundings suggests that the dense gas in post-shock stagnation regions is still the raw 
material out of which cores are built. We expect that in general core masses and collapse 
timescales will decrease with increasing density of the post-shock flow, which itself increases 
with increasing Mach number. Preliminary simulations of planar converging fiows that we 
have conducted indeed bear out this expectation, showing M^ore oc M~^. For planar con- 
verging fiows, many cores simultaneously grow and then collapse in the post-shock gas layer; 
unless this sheet were viewed exactly edge-on, the density jump at the shock front would not 
be apparent, and cores would be seen as surrounded by moderate-density gas. 

We also expect outside-in collapse followed by inside-out infall to be a generic feature of 
core evolution. Although the duration of this pressure-mediated collapse is ~ i//, it is unlike 
free-fall collapse in a crucial way: the core does not remain nearly-uniform. We speculate 
that the development of stratification during inside-out collapse will suppress growth of 
perturbations and sub-fragmentation for the non-spherical case. Even though the increase 
in density implies that the local Jeans mass becomes smaller and smaller, this is only true 
in the very center of the core. Instead, we expect that collapse of cores built within shocked 
converging fiows will produce single systems, which may be binary (or multiple) if the angular 
momentum is sufficient. Non-spherical converging fiows that create sheets and filaments of 
shocked gas and produce many such cores simultaneously could be the progenitors of stellar 
clusters. 

This research was supported by grant NNG-05GG43G from NASA. We are grateful to 
Phil Myers, Shantanu Basu, and the referee for careful reading and helpful comments on the 
manuscript. 
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Fig. 1. — Dashed curve, left scale: The temporal evolution of the mass accretion rate at 
inner edge of the grid, for collapse initiated from a critical BE sphere (see eq. Elfor units). 
Solid curve, right scale: Evolution of the central point mass (see eq. [TUfor units). Time is 
shown in units scaled by the density at the outer edge (see eq. [T2|) . 
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Fig. 2. — Density and velocity profiles {dotted curves) for converging- flow model with — 
1.05, at times r as noted in the upper panel for each pair. The radius is normalized by the 
central density, i.e. the abscissa is r^ATcGpcY^"^ /a. For the three upper panels on the right, 
the solid line is a fit to a BE sphere with temperature Tbe, with the fitted temperature noted 
in each panel. Dashed vertical lines denote the critical radius of a BE sphere with the same 
central density and sound speed a. Dot-dashed vertical lines mark twice this critical radius. 
The time for the leftmost pair is half of the collapse time Tcou- The time for the second and 
third pairs are when the shock reaches the critical BE radius and twice that value. The time 
for the last pair is the instant of collapse Tcou (defined in the simulations as Pc/po = 4 x 10^). 
The top-right panel shows with a dashed diagonal line the LP density profile D — 8.86^"^. 
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Fig. 3.- 



Same as Figure [2], for inflow Mach number Ma = 2.0. 
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Fig. 5.— 



Same as Figure [2], for inflow Mach number Ma = 7.0. 
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Fig. 6. — The shock front position versus time for Ma — 4. The sohd hne is the best fit to 
the hnear part, i.e. a constant-speed outward- propagating shock at early times. Collapse 
occurs for this model at Tcou — 0.16. 
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Fig. 7. — The speed of shock front at early times (asterisks), the radius of the core at 
the instant of collapse (diamonds), and the time at which collapse occurs (plus signs) as a 
function of Mq. Triangles show the ratio of icoul {ushTcoii), which is nearly constant, ranging 
from 0.34 to 0.42. The dotted line is the analytic estimate for Ush given in equation ( |23l) . 
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Fig. 8. — The density profiles (top panel) and the velocity profiles (bottom panel) at the 
instant of collapse for simulations with different Mach number Ma (as labeled). The thick 
dashed line in the top panel is the LP density profile D = 8.86.^"^. For higher Ma cases, the 
shock front is at smaller radius (in units of tq) and the post-shock speeds in the dense core 
are lower in magnitude. Units of length and velocity are given by equations f|T3|) and ffTTj) . 
respectively. Density is in units of the GMC ambient value. 
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Fig. 9. — The renormalized radius (top panel) and collapse time (bottom panel) of cores 
versus Ma- Diamonds show quantities normalized using the mean core density and aster- 
isks show quantities normalized using the Mach number (see text). The dot-dashed line 
in the lower panel shows the free-fall time at density Pmean in units {'^'KGpmeanY^'^ i 

T//(p^ean/P0)^/' = Vr(3/8)V2. 
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Fig. 10. — The density and velocity profiles during the accretion phase for = 1.05 (left) 
and Ma — 4 (right) models. Density profiles (top) and the velocity profiles (bottom) are 
each shown at three different instants: solid lines show the instant of core collapse, dashed 
lines show the instant when the gravitational rarefaction wave arrives at the shock front, and 
dot-dashed lines show the profiles at a late accretion stage. Numbers in the figure show the 
corresponding time for each instant. In the upper panels, the LP density profile D = 8.86^~^ 
is plotted for reference with thick dashed lines. The transitions from D oc (early) to 
D oc C"^/^ (late) and u (X (early) to m oc ^ (late) are evident. 
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Fig. 11. — The four stages of core evolution in the ideahzed models of this paper. (1) Core 
building via supersonic converging flow, yielding a stagnant, shock-bounded dense region. 
(2) Core collapse, propagating from outside to inside, leading to a density proflle p oc r^^. (3) 
Envelope infall, propagating from inside to outside, and resulting in free-fall onto the proto- 
star. (4) Late accretion of ambient gas onto the protostar. For real systems, evolution would 
be modified in several ways: converging flows would be non-spherical, angular momentum 
would lead to disk formation, and outflows would contribute to clearing the envelope. 
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Fig. 12. — Top panel shows the ratio of the duration of the supercritical collapsing stage 
ATsupcrit (asterisks) and the post-collapse infall stage ATinf (diamonds) to the prestellar 
lifetime Tcou of cores, as a function of Mach number. Bottom panel shows ratios Atsupcrit/tff 
(diamonds) and Atinf/tff (asterisks), which range from [0.8, 1], as a function of Mach 
number; here tff is computed using the mean density inside the shock at Tcou- 
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Fig. 13. — The temporal evolution of the mass accretion rate (thin solid lines; left axis) and 
the integrated mass of protostar (thick solid line; right axis) for models with Mach number 
Ma = 1.05, 2, 4 and 7, as labeled. Corresponding dashed lines show the results for models in 
which the inflow to the grid is suppressed after the end of the infall stage. The point marked 
i represents the end of the infall stage, p represents the instant when the density profile 
reaches p oc r^^/^ everywhere, q represents the instant when the velocity profile reaches 
V (X r~^/^ everywhere. Units of time, mass, and accretion rate are given by equations (fT2l) . 
( Tn|) and ([3]), respectively. 
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Fig. 14. — The protostellar mass at the end of the infall stage (diamonds), the core mass 
inside the shock at the instant Tcou of protostar formation (asterisks), and the critical BE 
sphere mass based on the post-shock density at time Tcou, all as a function of Mach number 
Ma- The mass unit is given by equation f|T^ . 



